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ABSTRACT 


The effect of solar or stellar radiation on dust particles’ trajectories (the 
Poynting-Robertson drag) has been studied by a number of authors and applied 
to interplanetary dust dynamics in numerical computations. Meanwhile some 
important features of dust flows can be studied analytically by implementing our 
novel hydrodynamical approach to use the continuity equation written in the 
particle’s orbital elements as coordinates (Gor’kavyi, Ozernoy, & Mather 1997). 
By employing this approach and integrating the continuity equation, we are able 
to hnd two integrals of motion when the Poynting-Robertson drag dominates 
the dissipative forces in the dust flow. These integrals of motion enable us 
to explore basic characteristics of dust flows from any sources in the Solar 
system (such as asteroids, comets, Kuiper belt, etc.) or in another planetary 
system. In particular, we have reproduced the classical solution n{r) cx r~^ that 
approximately represents the overall distribution of dust in the Solar system. 
We have also investigated various factors that could be responsible for the 
deviations of the power law index in n(r) oc from 6 = —1, including the 
influences of the orbital characteristics of dust sources, the evaporation of dust 
particles, as well as mixtures of dust particles of both asteroidal and cometary 
origin. We have calculated the masses and number densities of asteroidal and 
cometary components of the zodiacal cloud at different distances from the Sun. 


Subject headings: interplanetary medium - infrared: solar system 
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1. INTRODUCTION 


The Earth orbit is immersed in a multicomponent cloud of interplanetary dust and 
particles of various sizes. The structure and dynamics of this cloud are very important both 
for space missions and interpretations of numerous astronomical data. Small dust particle 
scattering (the zodiacal light) and emission contribute substantially to infrared data at 
AA 1 — 100 yum (Hauser 1995). 

The dynamics and evolution of interplanetary particles are determined by several 
basic effects which include: (i) the Poynting-Robertson (P-R) drag; (ii) resonance effects 
associated chiefly with Jupiter, Mars, Earth, and Venus; (iii) gravitational encounters with 
these planets, which occur in the form of elastic gravitational scattering of the particles by 
the planets; (iv) mutual collisions of the particles, and (v) evaporation of dust particles. As 
we have shown earlier (Gor’kavyi, Ozernoy, & Mather 1996, 1997), the dust flow evolution 
governed by those processes can be conveniently described by the continuity equation 
written in the space of orbital coordinates: 

dn{x^,t) ^ ^ ^ N^{xi,t) - N~{xi,t), (1) 

where i = 1...6, and Vi is velocity along the i-th axis in the phase space. In this equation, 
the div-term d{nvi)/dxi describes slow processes of change of particle orbital elements such 
as the particle transport due to both the P-R drag and the resonance effects. The terms 
N^{xi, t) and N~{xi, t) are responsible for fast processes such as the gravitational scattering 
of particles by the planets and the contact impacts of particles with each other and with 
the planets. 

Let us address the evolution of interplanetary particles governed by only two processes: 
(i) the P-R drag that continuously changes the particle’s orbit; and (ii) the gravitational 
scattering of particles by the planets, which changes the particle’s orbit like a jump (for 
detail, see Gor’kavyi, Ozernoy, & Mather 1997 referred to hereinafter as GOM). We 
introduce a small parameter e = Susc/dupR 1, where Supp and dusc are the changes, per 
unit time, due to the appropriate processes in the number density of dust particles at the 
typical point of the (a, e)-space, a being semimajor axis and e being eccentricity. This small 
parameter makes it possible to derive, by the Ghapman-Enskog approach, a set of equations 
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to get the first two consecutive approximations for particle density of the flow: 


dnpR{xi,t) 

dt 


d , 

+ ^ [npRVi) = 0 , 


dnsc{xi,t) 

dt 


N^{npR,Xi,t) - N {npR,Xi,t). 


( 2 ) 

( 3 ) 


Eq. (2) gives a first approximation to particle density associated with P-R drag, which is 
then used to get, with help of Eq. (3), a correcting term associated with the (second order) 
process of gravitational scattering. Eq. (3) applied to the dynamics of particles governed by 
gravitational scattering was analyzed in our previous paper (GOM). In the present paper, 
we emphasize Eq. (2) by exploring the Poynting-Robertson drag as the leading effect for 
the dynamics of small particles. The paper is organized as follows: In Ch. 1, we derive the 
integrals of motion of dust particles, assuming the PR-drag as the dominating effect. In 
Ch. 2, we hnd solutions of the continuity equation for the dust flow from a point source of 
dust (such as a comet or an asteroid). We also consider here the spatial distribution (in 
the coordinate space) of dust particles from such a source. In Ch. 3, we generalize those 
analytical solutions for evaporating dust particles, whose radii decrease with time. Ch. 4 
deals with dust flows produced by numerous real sources (217 comets and 5000 asteroids). 
In Ch. 5, we summarize our results and discuss the underlying assumptions as well as 
further work. 


2. INTEGRALS OF MOTION 


The rates Vi (i.e. da/dt,de/dt ,...) for the PR-drag have been given in a number of 
papers {e. g. , Burns et al. 1979; Lion, Zook & Jackson 1995): 

2 + 3e2 


da 

dt 

de 

dt 


= a- 


a (1 — e 
5e 


2C/2' 


= a 




( 4 ) 

( 5 ) 


where a = —M is the mass of the Sun, and (3 is the ratio of the solar radiation 
pressure force to the gravitation force: 


3LQ 


IQttG M c p I 


( 6 ) 
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where L is the solar luminosity, G is the gravitational constant, c is the speed of light, 
and I , p are the particle’s radius and density, accordingly. The coefficient Q depends on 
the optical properties of the particle and the ratio of its size to the wavelength of light; 
for a perfect particle Q = 1. For an asteroidal dust particle, we adopt 21 = 12 pm and 
p = 2.7 g/cm^, so that /? = 0.037 (e.g., Dermott et al. 1994). 

Along the particle’s trajectory, one gets 

<h ^ 5e(l - e^) 

da 2a(2 + 3e2)‘ ^ ^ 

By integrating this equation, the particle’s trajectory in a, e-space is easily obtained to be 

^ = const = Cl, (8) 

in agreement with Wyatt & Whipple (1950). This integral of motion implies a continuity of 
the particle’s trajectory in (a, e)-space. The trajectory is shown in Fig. 1. 


The flow of particles along evolutionary trajectories shown in Fig.l is governed by an 
equation (derived in the Appendix) that depends only on the variables a and e {i = const): 

f)ri f)ri 

ffla,e)— + Me)— = Mn,e), (9) 

where fl = 2a{2 + 3e^), /2 = 5e (1 — e^), and fz = n (6e^ — 1). The value of /s is given here 
for the case when the particle size does not depend on the orbital parameters a and e (for 
a more general case, see Appendix). The characteristics of the partial differential equation 
(9) (the integrals of motion) are Ci found above [Eq.(8)] and C 2 given by 

C 2 = ne^'^^Vl — e^. (10) 


The number density n of dust particles along the particle trajectory in a, e-space is shown 
as a function of semimajor axis a in Fig. 2 for various (7*2. The physical meaning of the 
integrals Ci and C 2 is the conservation of the particle’s flux along the flow under the P-R 
effect (note that the flux of the particles is conserved in conditions when the energy is not). 
Cl and C 2 as the integrals of motion for the particle’s flow under the P-R drag qualitatively 
play the same role as the Tisserand criterion does in the process of gravitational scattering 
(see GOM). As is shown in GOM, gravitational scattering results in “jumps” in the 
particle’s coordinates in a, e-space. However, those jumps obey the Tisserand criterion 
implying the conservation of the particle’s energy in the planet’s coordinate system under 
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the gravitational scattering, which is an elastic process. For a particle moving under 
the P-R drag alone, the energy is not conserved and the Tisserand criterion obviously 
breaks down, but the role of the integrals of motion is taken by Ci and 6 * 2 • Meanwhile 
for the particles which experience, along with the P-R drift, the (in fact, instantaneous) 
gravitational scattering by a planet, the Tisserand criterion holds as well. Indeed, during 
a short time interval of a gravitational scattering, the PR-drag practically does not change 
the particle’s energy. 

The integrals Ci ( 8 ) and C 2 (10) fully describe a stationary flow of dust particles 
produced by an arbitrary source of dust. As for the dust distribution in the coordinate 
space, it can be transformed from that in (a, e)-space with the help of Ci and C 2 by using 
the Hang (1958) double integral. This numerical procedure is rather complicated but, if 
we constrain ourselves to point sources in (a, e)-space (such as comets and asteroids with 
known a and e), the dust flows from such sources can be analyzed by means of a much more 
simple analytical technique to be presented in the next Section. 

3. PARTICLE NUMBER DENSITY DISTRIBUTION 
(Analytical Solutions for a Point Source) 


For a point source in (a, e)-space, two simple, 1-D continuity equations apply: 

d 

— (n{a)-Va) =0, 


da 

de 


( 11 ) 

( 12 ) 


;(n(e) ■ Z'e] ^ 0 . 

In orbital coordinates, 1-D distributions of the particle number density along the a-axis and 
e-axis are given, in the stationary state, accordingly by 


n{a) = C- 


a(2 + 3e2) ’ 


(13) 


Here C is a constant to be found from boundary conditions, and e in equation for n{a) can 
be eliminated with the help of the integral of motion Ci that describes, for a point source, 
the trajectory of motion of the entire dust flow. 
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A transform from the space of orbital elements (a, e) to r-space can be done using 
known analytical expressions (Sykes 1990). In r-space, the particles having orbital elements 
{a,e,f}, form a rather sophisticated torus-like cloud, whose number density is described by 
(Sykes 1990): 


P(r,0) = S{r) ■ 0(0), 


where 


S{r) = 




a?r , /g 2 _ ("I _ ’ 


0 ( 0 ) = 


(15) 

(16) 
(17) 


27r2v^cos^0 — cos^i 

with the limits a(l — e) < r < a(l -f e) and —i < 0 < i . For a volume density near the 
ecliptic plane, 0 = d and P(r, 0) oc (sin i)~^ . It is convenient to introduce the surface 
density 


ho 

(j{r) — J P{r, (p)dh, 


(18) 


-ho 


where dh = rd<p and therefore 


with 


a(r) = Co- ■ S'(r) ■ r 


I 

^ y 0(0)d0. 


(19) 


( 20 ) 


It is worth mentioning that functions (16) and (18) diverge at the pericenter and apocenter 
points. The total mass of dust particles having the same orbital parameters (a, e, i) in a 
thin spherical region between the radii (ri,r 2 ) is given by: 


12 

Am — 7r(ri + r2} J cr(r) dr. 


( 21 ) 


ri 


As a result of integration, one gets: 


A m(a, e) = 


vr(ri-Fr2)r . 1^2 


1 ,ri 


arcsin -(-1) — arcsin -(-1) 

e a e a 


( 22 ) 
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where 02(1 — 62 ) < r 2 ; ai(l + Ci) > ri . If the region under consideration, (ri,r 2 ), includes 
the pericenter or apocenter of the dust particle orbit, the integration should be done as far 
as the pericenter or apocenter, accordingly. Equation (22), as opposed to eqs. (16) or (18), 
does not involve divergences and is much more convenient in use. 

The distribution of the integrated surface density in the dust particle torus with 
different (a, e) is shown on Fig. 3. For a dust plume produced by any given source of dust, 
one can End, by using the solution (14) for n(e) of the continuity equation ( 12 ), the total 
dust mass in any region (ri,r 2 ) of the plume (taking into account the comment on the 
integration limits) is: 

62 

A M = J n{e) ■ Am{a, e) de. (23) 

61 

The result is given by 

AM = Cm- 7r(ri + r2) / ^^ 

ei 

f 1 rr2(l — e^) ] lrri(l —e^) 1 1 

I eL Cie^A J Cie^A Jj ^ ^ 

The integration limits are easily found by comparing the location of the apocenter or 
pericenter of particle orbit relative to the boundaries of the region under consideration, 
(ri,r 2 ). Thus the problem of a transform from the space of orbital elements to the r-space 
is solved by reducing it to a comparatively simple single integration to be done numerically. 

The above technique can be applied to both stationary and non-stationary flows. The 
source of dust can appear owing to an instantaneous ejection, e.g. be associated with a 
collision of asteroids (Sykes 1990; Mann et ah 1996). The resulting drift of dust particles 
(assumed to have the same size and located initially as a point source at Oq, Cq) is shown, 
for three instants of time, in Fig. 3. 

Let us address a situation when the source of dust has recently appeared and continues 
to operate. This happens, for example, when a new comet appears in a given point of 
(a, e)-space as a result of a ‘jump’ due to the strong gravitational scattering of that comet 
by a planet. Suppose that the first dust particles ejected from such a source have reached 
the semimajor axis Omm- In order to derive the density distribution in the formed dust 
cloud, it is sufficient to integrate equations(23) and (24) down to Cmin, and not through 
the entire dust plume. The value of e^m can be found from eq. ( 8 ) at the known amin- 
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The appearance of this dust plume is shown, as a function of time, in Fig. 4. Fig. 5 
demonstrates a steady-state plume formed from various sources of dust particles. As can 
be seen, small eccentricities result in the surface density a{r) ^ const, which implies that 
the volume density (i.e. accounting for inclination i) is n(r) oc r~^. This result reproduces 
the well-known characteristic density run of the interplanetary dust. If the source of dust 
has an eccentric orbit (such as the Encke comet), the density run is much steeper. The 
boundary between the plateau and cut-off of all three curves on Fig. 5 is located at the 
pericenter of the source of dust, q = a(l — e), which is due to the fact that the pericenter of 
a particle experiencing the P-R drag undergoes the least shift compared to its semimajor 
axis and apocenter. 

4. INTEGRALS OF MOTION FOR EVAPORATING PARTICLES 


Here we deal with a situation, important in the dynamics of the interplanetary 
particles, when one needs to take into account an evaporative effect of the solar wind upon 
a particle (see e.g. Johnson 1991). We describe the mass loss of a particle of mass m by 

^ = (25) 


where m* = fXsmp is the mass of molecules sputtered from the dust particle’s surface, fUg is 
the molecular weight, rup is the proton mass, Ni is the number of particles of the Tth kind 
in the solar wind colliding with the dust particle, and Yi is the sputtering coefficient, i.e. 
the number of moleculs sputtered from the dust particle by one wind particle. If the dust 
particles that drift to the Sun owing to the PR-drag have small eccentricities, the particle’s 
size, /, changes according to a simple law I oc a^, where the power law index A is given by 
(Shestakova 1994): 


A 


Jq l^s Y 
6(1 -|- Jo) Hw 


(26) 


Here is an average molecular weight of the solar wind, Jq = Mc^/L = 0.3 is the ratio of 
the corpuscular luminosity of the Sun to the radiative one L, and it is assumed that Q = 1. 
If one only accounts for the contribution of protons and a-particles to the dust particle’s 
mass changes, one hnds (Shestakova 1994): 


A = 0.032/isYs, 


(27) 
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which yields A = 0.025 for silicates, and A = 0.56 for ice (H 2 O). Bearing in mind that both 
f3 and a oc 1// oc a“^, equation (9) takes the form (see Appendix): 

(IT! (IT! r ~\ 

2 a(2 + 36^)— + 5e(l - e^)— = n 66^(1 + A) + A - 1 . (28) 

oa ae L J 

From this equation, the first integral given by eq. (8) is left without any changes, which 
implies that the dust particle subjected to evaporation continues to move along the same 
a, e-trajectory. Meanwhile the second integral given by eq. (10) is modified into the 
following equation (see the Appendix): 

C'2 = nei/'^-("/5)"(l-e2)"+V2_ (29) 

Such a modification results from an increase of the particle’s drift velocity in the (a, e)-space 
owing to a decrease of the particle’s size due to evaporation. 

The effect of the particle’s evaporation can also be accounted for while deriving the 
densities n(a) and n(e) in eqs. (13) and (14). In order to determine the mass of dust in 
any region (ri,r 2 ), one needs to account for the dependence of a upon the function l{a,e) 
in the integrand of eq. (24). Besides, one needs to account for a decrease in mass of the 
evaporating dust particle. As a result, a factor (or a^^), with a corresponding change 
in the constant Cm, appears in the integrand of eq. (24). Fig. 6 demonstrates a ‘plume’ 
formed from a typical asteroid for different A = 0; 0.05; 0.125; 0.25. 

5. SOURCES THAT FEED THE ZODIACAL CLOUD 

As the sources of dust, we have incorporated the first 5000 numbered asteroids (ITA 
1994), which comprise several distinct families (Zapalla et al. 1995) as well as 217 comets 
with e < 1 taken from Marsden & Williams (1992). The positions of these sources on the 
diagram (sinqg) are presented on Figs. 7 (asteroids) and 8 (comets). One can see that 
the asteroid percenters tend to concentrate in the two regions 1.7-2.15 AU and 2.6-2.9 
AU. The volume density n{r) = a/{rsini), where rsini is the thickness of the dust plume. 
Therefore the lower the position of an asteroidal group in Fig. 7 (small sini), the larger is 
its contribution to n{r) in the ecliptic plane. As for comets, their pericenters are distributed 
more or less uniformly at 0.5A5'A2.5 AU. 
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In order to compute the contributions of asteroidal and comet sources to the volume 
density of the interplanetary dust, we adopt the simplest assumptions: (i) all 5000 asteroids 
produce dust particles with the same rate and (ii) comets produce dust particles whenever 
they pass pericenter so that their contribution into the zodiacal cloud is oc l/(Tg'’'), where 
T is the comet orbital period and typically 7 = 2 . We neglect the comets (sungrazers) with 
a(l — e) < 0.01, assuming that the dust produced so close to the Sun evaporates very soon. 

To visualize easily any deviations of n(r) from the ‘standard’ solution n{r) = r~^ 
while making comparison with the data on the volume density distribution of 
interplanetary/zodiacal dust, it is convenient to deal with n{r)r = a/sini. The 
expected contribution of the above 5000 asteroidal sources to n{r)r is shown in Fig. 9. It 
can be seen that the cut-off at r ~ 1.7 — 2.0 AU, associated with the asteroidal distribution 
in q, is mainly due to the asteroids belonging to the inner and middle groups. We have also 
computed n(r) due to the 1000 largest asteroids with a result that differs insignihcantly 
from the entire sample. Therefore, our neglecting the asteroidal distribution in sizes is not 
critical for our conclusions. The expected contribution of the above 217 comet sources to 
n{r)r is shown in Fig. 10. It can be seen that the density run does not have kinks at r > 0.8 
AU and is much steeper compared to what is due to the asteroidal sources (see Fig. 9). 

Let us consider the combined cometary and asteroidal contributions into the 
interplanetary dust density. The total value of n{r)r for a mixture of 25% (in number 
density) cometary and 75% asteroidal components is shown in Fig. 11. 

The density run has been approximated by n(r) oc r^. The power law indices 6 at 
different r are given in Table 1 for a particular 7 = 2 . 

%% EDITOR: Put Table 1 here %% 

We have found that the results for 7 = 1 and 7 = 8 are qualitatively similar: the 
ratio of the volume densities (asteroids/comets) has been found to be (65%/35%) at 
7=1 and (80%/20%) at 7 = 3. Accounting for evaporation of dust particles changes the 
results in a more substantial way: if one takes Aa=0.025 for asteroidal dust and Ac=0.56 
for cometary dust (Shestakova 1994) then one gets the power law index 6 = —1.3 in the 
region 1.0 < r < 1.7 AU when the ratio of the volume densities (asteroids/comets) equals 
45%/55%. 

The asteroidal/cometary mass ratios, which provide the observed 6 = —1.3 for different 
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possible values of 7, A^, and Ac, are given in Table 2. 

%% EDITOR: Put Table 2 here %% 

The mass ratios presented in Table 2 for spherical layers at various r differ from the 
ratio of volume densities in the ecliptic plane (the 4th line of Table) because of different 
inclinations of the cometary and asteroidal sources of dust and, therefore, of different 
thickness of the dust cloud. The last column of Table 2 is fairly consistent with the 
estimations for asteroidal/cometary ratios obtained from the IRAS observations while 
modelling the shape of the zodiacal cloud (Liou, Dermott, & Xu 1995). It is worth noting 
that the evaporative power law index Aq influences those results insignihcantly: if one puts 
Aa=0.0, the data in the last column of Table 2 would change by no more than 1%. 

6. DISCUSSION AND CONCLUSIONS 


The continuity equation written in the space of orbital elements, which we have 
proposed to use as a tool for studying the dynamical evolution of interplanetary particles, 
is very effective not only for gravitational scattering (Gor’kavyi, Ozernoy & Mather 
1996, GOM) but, as is shown above, for the chief factor of dynamical evolution of the 
interplanetary dust cloud - the Poynting-Robert son effect. By solving the continuity 
equation, we have found two integrals of motion governed by the P-R drag and, as a result, 
have derived the density run for a quasi-stationary state of the interplanetary dust cloud. 

Based on our approach, we are able to compute the formation of a dust cloud from 
many cometary, asteroidal, etc. dust sources, either observed or assumed to exist, by 
adopting various reasonable assumptions about the rates of dust production. As an 
outcome, we get the theoretical model for dust distribution both in radius and latitude, in 
the Solar system. We call ‘the reference model’ the model that only accounts for the solar 
gravitational held and the P-R effect . By accounting for the other dynamical factors, we 
can get corrections to the reference model. 

Gomputations performed in the framework of the reference model yield radial density 
runs of the interplanetary dust obtained nnder various model assumptions and described 
above. They enable us to derive some general conclusions as follows: 


The asteroidal component of the interplanetary dust is typically characterized by 
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n{r) oc r~^ up to 1.7 AU with a sharp cut-off at r72 AU [e.g n{r) oc at r ~ 3 
AU]. The predicted cut-off is entirely consistent with what is observed for the zodiacal 
cloud and the main core-population of the dust in the Solar system (Divine 1993); 

• The power law index 5 in the distribution of the volume density of the zodiacal dust 
particles, n{r) oc r'^, depends on several factors: 

— at radii smaller than the pericenter of the source of dust, (r < a(l — e)), the 
index 5 = (5i differs signihcantly from that, 6 = S 2 , a.t radii larger than the 
pericenter of the source of dust, (r > a(l — e)); 

— both (5i and S 2 depend strongly on eccentricity of the source: as e increases from 
0 to 1, (5i increases from —1 to 0 whereas the value of S 2 increases from —00 to 
-2.5; 

— accounting for the particle’s evaporation, both (5i and 62 tend to increase (see 
Fig. 6); 

— for a system of dust sources distributed in a and/or e, the index S differs 
substantially inside and outside the system (see Fig. 12). For such distributed 
sources, there is a smooth transition between (and not a jump from) (5i and S 2 
near the pericenter. 

• In the region r < 1.7 AU, observations yield n{r) oc (Divine 1993) or 
(Levasseur-Regourd 1996). The observed deviation from the classical law r~^ may be 
related to the contribution of the cometary dust component. This component (which 
alone would result in n{r) oc r“^-^) is expected to contribute as much as about 25% in 
number density (or about 50% in mass) to the dust cloud near Earth; 

• The ratio of the cometary dust contribution to that of the asteroidal dust into n(r) 
has a minimum at r ~ 1.8-2 AU and it is large in the inner region of the Solar system 
(a < 1 AU) as well as in the outer region of the main asteroid belt (a > 3 AU). 

• The dust created by the Kuiper belt should £11 the inner region of the Solar system 
as a more or less homogeneous layer with a{r) = const. At r > 3 — 4 AU that dust 
may dominate over the comet and asteroidal dust components. 

Let us turn to discussing how the adopted approximations influence our conclusions. 
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The basic approximations of our reference model are: (i) a steady-state character of the 
model and (ii) neglect of other dynamical factors except the P-R drag. 

Approximation (i) could be easily abandoned - we are able to compute the non¬ 
stationary processes as well. For instance, Fig. 4 shows different stages of dust plume 
formation of a source that once appeared and continues to produce dust particles. Therefore, 
this approximation does not imply any restriction in our approach or a failure of our model; 
it only results from an incompleteness of our knowledge of whether the real sources of dust 
are stationary or not. 

Approximation (ii) implies that our solutions fail to describe such phenomena as the 
resonance circumsolar dust rings near Earth (Jackson & Zook 1989), which may contain as 
many as ~ 10% of the particles of the zodiacal cloud (Dermott et ah 1994; Liou, Dermott, 
& Xu 1995). In our further work, we will fully incorporate the resonance and related effects. 

Some minor approximations in our computations, such as the intensity of comets as 
dust sources as a function of their distances from the Sun, and the law of dust particle 
evaporation, are also worth rehning. Besides, we have ignored in the above computations 
the distribution of dust sources in size and chemical abundance as well as the effects of 
observational selection in asteroidal and cometary populations. However, all uncertainties 
in the those efffects can be easily explored by varying the parameters of the model for the 
zodiacal cloud formation. For instance, while computing a multi-component model of the 
dust cloud, we need to compute the reference model and get the structure of the dust cloud 
consisting of the particles of a given size. The value of f3 influences not just the speed of a 
particle’s motion along the orbital trajectories but, when (3 is not too small compared to 1, 
one needs to account for the change of orbital coordinates of the dust particles beginning 
with the instant of their escape from the source (see, e.g., Jackson & Zook 1992). 

Despite the above approximations, our reference model has many principal advantages. 
Being semi-analytical it is flexible enough to be able to compute the dust cloud formation 
from a large number of various sources having an arbitrary distribution in size, composition, 
etc. It takes a short CPU time (e.g., 4 min for 217 comets and 1 hour for 5000 asteroids 
with a 486/DX4-100) to compute the typical reference model. Therefore it is feasible to deal 
with multi-variant dust population sources while fitting the available cometary, asteroidal, 
and far infrared COBE data. Then one could employ our basic model of the dust cloud to 
incorporate the other effects such as gravitational scattering of dust particles by planets 
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(Gor’kavyi, Ozernoy & Mather 1997), resonance capture of dust (Jackson & Zook 1989; 
Dermott et al. 1994), particle multiple collisions, etc. 

We conclude by emphasizing that our continuity equation (1) is well suited to the study 
of various large scale properties of the interplanetary dust cloud. Our approach might be 
fruitful for solving a number of related problems, including the structure of protoplanetary 
and circumstellar disks. 
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Appendix 

Integrals of motion. 


Here we derive the integrals of motion based on eq. (11). At first, we address the 
case when the particle’s size does not depend upon the orbital parameters, and thus the 
coefficient a in eqs. (4) and (5) does not depend upon the semimajor axis a and eccentricity 
e. By differentiating eqs. (7) and (8), one hnds dva/da and dve/de, respectively: 

2 + 3e2 


dva 2 + 3e^ d /l\ 

^ “ “(1 _e2)3/2^1^aj “ 


.2^3/2^ 


dve 5 (9 / e 

^ ^ "^^1^(1 -e2)V2^ 


= a- 


2a2(l-e2)3/2' 

By substituting eqs. (7), (8), (Al), and (A2) into eq. (11), one gets: 


(.41) 

(A2) 


q;(2 + 3e^) dn na{2 + 3e^) 


+ 


Soe 


dn 


j(l_e2)3/2 5a 2a‘^{l - de 

5na 


2 a 2 (l 


. 3/2 


= 0 . 


Having multiplied eq. (A3) by 2q; ^a^(l — one hnds dually: 

2a(2 + + 5e(l - e^)^ = n{Qe^ - 1) 

oa oe 

This partial differential equation implies the following ordinary differential equation : 

da de dn 

2a{2 + 3e2) “ 5e(l - e^) “ n(6e2 - 1) 


(.43) 


(M) 


(.45) 


Each pair of equations comprising eq. (A5) is solved by separating the variables to give 

(A6) 


1 2 /• 2 + 3e2 

57 e(l-e2)'^'' 


and 


1 , 1 /• 6 e 2 - 1 , 


(.47) 


The integrals can be easily calculated to give the integrals of motion discussed in the main 
text: 


C'i = 


a(l — e^) 
eTs 


(A8) 
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and 

C 2 = (A9) 

Let us now address the case when the particle’s size, owing to evaporation, depend 
on the orbital parameters a and e. For simplicity, we assume that the particle’s radius I 
changes with a according to a power law. Then the coefficient a oc oc a~^, i.e. 

a = aQa~^. (^ 10 ) 


In this case, one gets, instead of eq. (Al): 


dvg 

da 


2 + 3e^ + 2 + 


(1 - da\a 
as well as the following equation, instead of eq. (A2): 


a2+A(l _ £2)3/2’ 


(All) 


dVe 

de 


San 


= a 


2a2(l - 2a^+\l - e^) 


3/2 ■ 


(Aii) 


Having substituted eqs. (AlO), (All), (A12), (7), and (8) into eq. (11), one gets after some 
algebra 

2 a(2 + 3e^)^ + 5e(l — e^)^ = n 6e^(l + A) + A — 1 . (A13) 

da ae L -i 

Since the l.h.s. of eq. (A13) and (A4) are equal to each other, Ci is kept the same as in eq. 
(A8) for a non-evaporative case. Meanwhile for the r.h.s. of eq. (A13) one gets, similar to 
eq. (A7): 

r ^ T r Xl -U X - 1 

(3114) 


1 1 /" 6e^(l + A) + A — 1 

^" = 5 /— 


which yields the second integral of motion for (evaporating) particles of a changing radius: 


(3115) 


C 2 = „el/5-(4/5)A(i _ e2)A+l/2^ 
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Table 1. Power Law Indices 6 


Region 

Asteroidal dust 

Cometary dust 

75% ast. dust ± 25% com. dust 

r = 0.5 - 1.7 AU 

-1.03 ±0.02 

-2.17 ±0.38 

-1.29 ±0.06 

r = 1.0 - 1.7 AU 

-1.04 ±0.02 

-2.36 ±0.29 

-1.29 ±0.03 

r ^ 2.0 AU 

-2.6 

-2.8 

-2.6 

r ^ 2.5 AU 

-4.0 

-2.9 

-3.8 

r ^ 3.0 AU 

-7.4 

-3.1 

-6.4 


Table 2. Mass of Asteroidal Dust / Mass of Cometary Dust (% / %) 


r(AU) 

1 

7 = 1 

Aa = 0 

Ae = 0 

Ua/ric = 65/35 

7 = 2 

Aa = 0 

Ac = 0 

na/ric = 75/25 

7 = 3 

Aa = 0 

Ae = 0 

UalUc = 80/20 

7 = 2 

Aa = 0.025 

Ae = 0.56 
na/ric = 45/55 

0.0 - 1.0 

39/61 

42/58 

37/63 

19/81 

1.0-3.0 

56/44 

69/31 

75/25 

28/72 

3.0-5.0 

21/79 

33/67 

40/60 

5/95 

0.0-5.0 

46/54 

54/46 

52/48 

21/79 

0.95-1.05 

46/54 

58/42 

62/38 

25/75 


*) asteroidal/cometary volume densities of dust at r = 1 AU in the ecliptic plane. 
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Figure Captions 

Fig. 1. Changes in the particle’s position in the (a, e)-space due to the Poynting- 
Robertson drag. Initial semimajor axis is Oq = 2.5 AU; initial eccentricities are 
Co = 0.15, 0.29, 0.43, 0.57, 0.71, 0.85, and 0.99 from bottom to top, accordingly (the 
corresponding Ci are 11.15, 6.16, 4.00, 2.64, 1.63, 0.79, and 0.05). The transforms of orbits 
A, B, and C to the r-space are shown in Fig. 3. 

Fig. 2. The number density of dust particles along the particle trajectory in 
a, e-space, n, is shown as a function of semimajor axis a for various initial eccentricities 
Co = 0.15, 0.29, 0.43, 0.57, 0.71, 0.85, and 0.99 (with no = 1). The upper dashed line, 
n = corresponds to the circular case (cq = 0 ); the lower dashed line, n = C4 
corresponds to the parabolic case (cq = 1 ). 

Fig. 3. The surface density, o'(r), integrated over small (r 2 — ri = 0.01 AU) regions 
as a function of radius r for three different points A, B, and C shown in Fig. 1. Initial 
(ao,eo)-data for those points are, respectively: (2.5, 0.29); (1.3, 0.14); and (0.6, 0.054). 

Fig. 4. Seven stages of the evolution for the surface density of a cometary plume 
having qq = 3.0, cq = 0.5. Curves labeled 1 to 7 corresponds to the instants of time by 
which the semimajor axis of particles decreases by jttQ. 

Fig. 5. The surface density a as a function of radius r for the a, e-curve of three different 
sources of dust having initial semimajor axis Oq = 2.22 and eccentricities Cq = 0.15; 0.5; 0.85 
(curves 1, 2, and 3, accordingly). The first case is typical for the inner asteroids, and the 
last case is the Encke comet. 

Fig. 6. The surface density a as a function of radius r for the a, e-curve of sources of 
dust with oo = 2.22 and eccentricity eo = 0.15 at different A = 0; 0.05; 0.125; and 0.25 
(curves 1, 2, 3, and 4, accordingly). 

Fig. 7. 5000 asteroids (taken from ITA 1994) plotted in the [sini, a(l — e)]-space. 
Enhanced concentrations of points indicate at least 4 to 5 different asteroidal groups. 
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Fig. 8. Comets in the [sini, a(l — e)]-space (Marsden & Williams 1992). Notations: 
crosses for 97 short-periodic comets with a < 20 AU and i < 90°; squares for 67 long-periodic 
comets with a > 20 AU and i < 90°; triangles for 53 retrograde comets with i > 90° (all but 
three comets have a < 20 AU). 

Fig. 9. Function n{r)r = u{r )/sini as a function of radius r for dust particles of 
asteroidal origin. Curves 1, 2, and 3 correspond, respectively, to 1816 asteroids of the inner 
zone (a < 2.501 AU); to 1271 asteroids of the mid zone (2.501 < a < 2.825 AU); and to 1913 
asteroids of the outer zone (a > 2.825 AU). Curve 4 is the integral for all 5000 asteroids. 
One can see that n{r)r has a cut-off at r ~ 1.7 — 2.0 AU. The curve consisting of squares 
corresponds to the dust produced by all 5000 asteroids with accounting for evaporation 
(A = 0.025). 

Fig. 10. Function n{r)r = u{r )/sini as a function of radius r for dust particles of 
cometary origin. Curves 1, 2, and 3 correspond, respectively, to all comets with 7 = 1, 2, 
and 3. The solid curve (the evaporative case) corresponds to the dust produced by comets 
with 7 = 2 and A = 0.56. 

Fig. 11. Functions n{r)r for cometary dust (squares) and asteroidal dust (triangles). 
Each function is normalized to unity at r = 1 AU. Crosses stand for a mixture of cometary 
component (25%) and asteroidal component (75%) (%% are taken in number density at 
r = 1 AU in the ecliptic plane). 

Fig. 12. Functions S{r) for a single comet with a = 2.22AU and e = 0.85 (curve 1); 50 
hctitious comets uniformly distributed in e within the range of 0.6 < e < 0.85 (curve 2); a 
sample of 217 real comets (curve 3); a single asteroid having a = 2AU and e = 0.15 (curve 

4) ; 200 hctitious asteroids uniformly distributed in a within the range of 2 < a < 3 (curve 

5) ; and a sample of 5000 real asteroids (curve 6). 
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